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Abstract 

An  essential  need  of  the  US  Air  Force  is  the 
discovery ,  development ,  and  fielding  of  new ,  energetic 
materials  for  advanced  chemical  propulsion  in  space  and 
missile  applications.  Some  of  the  key  factors  driving  the 
requirement  for  new  chemical  propellants  include:  (a) 
improved  performance  in  terms  of  increased  specific 
impulse  and  density ,  (b)  reduced  sensitivity  to  external 
stimuli  such  as  impact,  friction,  shock,  and  electrostatic 
discharge,  and  (c)  mitigation  of  environmental  and 
toxicological  hazards  (and  the  resulting  costs )  associated 
with  currently  used  propellants. 

A  class  of  compounds  which  can  potentially  meet 
these  requirements  is  known  as  ionic  Liquids  (ILs),  which 
are  chemical  salts  with  unusually  low  melting  points.  The 
physical  and  chemical  properties  of  ILs  render  them 
useful  for  many  purposes,  most  notably  as 
environmentally  benign  ((< green  ”)  solvents/reaction 
media  but  also  as  catalysts,  electrolytes,  etc.  From  a 
Department  of  Defense  (DoD)  perspective,  ILs  are  being 
explored  as  new  propellants,  explosives,  and  munitions. 
The  Air  Force,  in  particular,  is  interested  in  ILs  as 
potential  replacements  for  currently  used 
monopropellants  such  as  hydrazine  -  which  is 
carcinogenic,  highly  toxic,  and  has  relatively  modest 
performance  characteristics.  In  contrast,  many  ILs  have 
superior  densities  and  specific  impulses  as  well  as 
significantly  reduced  sensitivity  and  toxicity 
characteristics.  Furthermore,  their  properties  can  be 
carefully  tuned  via  the  choice  of  the  component  ions. 

The  overall  objective  of  the  Design  of  Energetic  Ionic 
Liquids  challenge  project  is  to  address  several  key 
technical  issues  and  challenges  associated  with  the 
characterization,  design,  and  development  of  ILs  as  new 
monopropellants.  Among  these,  for  example,  are  a 
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fundamental  understanding  of  the  (in)stability  of  ILs,  the 
intrinsic  nature  of  the  short-  and  long-range  structure  and 
interactions  between  the  component  ions,  and 
identification  of  the  key  steps  in  the  initial  stages  of 
decomposition  and  combustion.  A  hierarchy  of 
computational  approaches  is  employed,  including 
atomistic,  high-level  quantum  chemical  methods  applied 
to  individual  ions  and  ion  clusters,  condensed  phase 
atomistic  molecular  dynamics  simulations  utilizing 
polarizable  force  fields,  and  me  so  scale -lev  el  simulations 
of  bulk  ionic  liquids  based  upon  multiscale  coarse 
graining  techniques. 

1.  Introduction 

The  design  of  new  high  energy  density  materials, 
which  are  more  efficient,  reliable,  and  environmentally 
benign  than  existing  rocket  propellants,  is  a  high  DoD 
priority.  The  focus  of  this  effort  has  been  on  the 
development  of  new  propellants  and  energetic  additives, 
including  highly  strained  hydrocarbons,  poly  nitrogen 
compounds,  and  advanced  monopropellants.  Some  of  the 
issues  that  must  be  addressed  in  theoretical  efforts  to 
design  new  energetic  materials  include  an  assessment  of 
their  energy  content,  their  thermodynamic  and  kinetic 
stability,  and  the  design  of  new  synthetic  routes  to 
proposed  new  compounds  that  have  not  yet  been 
synthesized. 

A  specific  area  of  interest  to  the  DoD  is  the  discovery 
of  a  suitable  replacement  for  hydrazine,  a  widely  used 
monopropellant  for  low-thrust  propulsion  applications 
such  as  orbital  maneuvering  and  satellite  stationkeeping. 
The  desire  to  replace  hydrazine  is  motivated  by  several 
factors.  Perhaps  the  most  severe  limitation  of  hydrazine 
is  its  carcinogenic  nature  and  extreme  respiratory  and 
dermatological  toxicity,  with  correspondingly  large  costs 
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associated  with  controlling  these  environmental  and 
biological  hazards.  Furthermore,  the  performance  of 
hydrazine  as  a  monopropellant  is  rather  modest  due  to  its 
relatively  low  density  and  specific  impulse  compared  to  a 
prototypical  ionic  monopropellant  salt  such  as  4 -amino- 
1, 2, 4-triazolium  dinitramide.  The  replacement  of 
hydrazine  with  more  energetic,  less  hazardous  energetic 
monopropellants  is  clearly  needed. 

A  specific  type  of  energetic  material  of  current 
interest  is  derived  from  a  broad  class  of  compounds 
known  generically  as  ionic  liquids  (ILs),  which  are 
chemical  salts  with  unusually  low  melting  points;  e.g., 
below  100°  C.  The  general  interest  in  ILs  has  focused 
mainly  on  their  use  as  environmentally  benign  (“green”) 
solvents  for  a  wide  range  of  chemical  reactions.  Some  of 
the  properties  of  ILs  which  make  them  attractive  as 
solvents  include  their  low  vapor  pressure,  large  liquid 
ranges,  and  thermal  stability.  The  interest  in  ILs  as  new 
monopropellants  stems  from  several  factors.  For 
example,  the  properties  of  ILs,  including  their  energy 
content,  can  be  “tuned”  through  a  judicious  choice  of 
component  ions  and  their  substituents.  Furthermore,  the 
virtually  nonexistent  vapor  pressure  of  ILs  greatly  reduces 
the  environmental  and  biological  hazards  due  to 
respiratory  and  dermatological  toxicity.  Finally,  the 
densities  of  ILs  generally  are  significantly  greater  than 
those  of  conventional  liquid  monopropellants  such  as 
hydrazine. 

Although  there  have  been  extensive  experimental 
studies  of  chemical  reactions  in  ILs,  little  has  been  done 
in  the  area  of  characterization  of  the  fundamental 
chemical  and  physical  properties  of  ILs.  In  particular, 
one  of  the  most  pressing  needs  in  the  broader  area  of  IL 
development,  and  particularly  in  the  design  of  energetic 
ILs,  is  the  application  of  robust  theoretical  methods  for 
the  reliable  prediction  of  IL  heats  of  formation,  synthesis 
routes,  phase  transitions,  ion  conformations,  thermal 
stabilities,  densities,  and  viscosities.  The  focus  of  this 
study  is  on  the  characterization,  design,  and  synthesis  of 
the  next  generation  of  monopropellants  for  rocket 
propulsion  applications. 

2.  Computational  Methods 

An  integrated  approach  utilizing  multiple 
computational  methods  is  used  to  predict  and  characterize 
the  intrinsic  and  bulk  properties  of  energetic  ionic  liquids. 
At  the  molecular  level,  highly  accurate  electronic 
structure  methods  are  used  to  predict  the  fundamental 
properties  of  the  ionic  liquid  components,  including 
molecular  structures,  charge  delocalization,  heats  of 
formation,  and  proton  transfer  reaction  pathways  and 
barriers.  Geometries,  electronic  structures,  and  properties 
(including  heats  of  formation)  of  the  component  ions  are 


predicted  using  second-order  perturbation  theory  (MP2, 
also  known  as  MBPT(2)[1]),  density  functional  theory 
(DFT)[2],  coupled  cluster  theory  (CCSD(T)[3])  and  the 
“Gaussian-N”  (GN)[4]  methods.  The  Nuclear-Electronic 
Orbital  (NEO)[5]  approach  is  used  for  capturing  the 
quantum  dynamical  effects  of  hydrogen  bonding  and 
proton  transfer.  In  addition,  the  complex  spectrum  of 
ionic  liquid  physical  properties  requires  utilization  of 
polarizable  force  fields  in  condensed  phase  simulations  in 
order  to  obtain  reliable  predictions  of  many  key  bulk 
properties. 

3.  Results  and  Discussion 

Cation-cation  n-n  Stacking'.  A  theoeretical  study[6]  of  the 
gas-phase  geometries  and  interactions  in  small  clusters  of 
ions  suggests  the  presence  of  n-n  stacking  structures 
between  the  cations  in  dimers  of  1,2, 4-triazolium 
dinitramide  (two  cations  and  two  anions)  as  well  as  1,2,4- 
triazolium  chloride.  These  types  of  structures  can  be 
formed  in  dielectric  solvents  due  to  screening  of  the 
charge-charge  repulsion  forces, [7]  but  such  species  in 
small  structures  have  not  been  previously  reported. 

Of  the  numerous  structures  found  for  the  two  pairs  of 
1,2, 4-triazolium  and  dinitramide  ions,  or  the  pairs  of 
corresponding  neutral  1,2,4-triazole  and  dinitramine 
molecules,  the  most  stable  MP2/aug-cc-pvdz[8]  optimized 
geometries  are  shown  in  Figure  1.  In  the  ionic  structure, 
each  1,2, 4-triazolium  forms  two  hydrogen  bonds,  via  the 
hydrogens  on  the  N  atoms,  to  the  O  atoms  of  the 
dinitramide  ions.  Interestingly,  this  structure  exhibits 
parallel  stacking  of  the  two  cationic  1,2, 4-triazolium 
rings.  The  interplane  distance  is  -3.2  A,  with  a  parallel 
displacement  of  ~1.4  A.  The  corresponding  neutral 
tetramer  shows  a  similar  parallel  stacking  arrangement  of 
the  triazole  rings. 

Furthermore,  it  is  of  interest  to  determine  the  cluster 
size  at  which  the  ion  pair  structures  become  more  stable 
than  the  corresponding  neutral  pair  structures.  A  previous 
study  predicted  that  ion  pair  dimers  are  typically  higher  in 
energy  than  neutral  pair  dimers. [9]  Including  zero  point 
vibrational  energy  (ZPVE)  corrections,  the  ionic  tetramer 
in  Figure  1  is  1.2  kilocalorie  /  mole  (kcal/mol)  lower  than 
that  of  the  neutral  one.  The  MP2  method  tends  to  predict 
higher  energies  for  ionic  species  vs.  neutral  species, [9]  so 
more  accurate  CCSD(T)/aug-cc-pVDZ  energy 
calculations  of  these  two  tetramer  structures  were  desired. 
However,  since  the  computational  cost  of  CCSD(T)/aug- 
cc-pVDZ  is  prohibitive,  these  energies  were 
approximated  from  the  MP2/aug-cc-pVDZ  energies  by 
estimating  the  electron  correlation  energy  differences 
using  three  independent  methods:  (1)  the  differences 
between  the  MP2/cc-pVDZ[8]  and  CCSD(T)/cc-pVDZ 
energies  of  the  tetramers,  (2)  the  differences  between  the 
MP2/aug-cc-pVDZ  and  CCSD(T)/aug-cc-pVDZ  energies 
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of  the  twelve  pairs  of  dimers  in  these  two  tetramers,  and 
(3)  the  differences  between  the  MP2/aug-cc-pVDZ  and 
CCSD(T)/aug-cc-pVDZ  energies  of  the  eight  monomers 
in  these  two  tetramers.  Using  these  three  methods,  and 
including  ZPVE  corrections,  the  estimated  CCSD(T)/aug- 
cc-pVDZ  energy  of  the  ionic  tetramer  is  lower  than  that  of 
the  neutral  tetramer  by  5.7,  7.3,  and  7.7  kcal/mol, 
respectively. 

The  driving  force  for  the  formation  of  the  unusual 
cation-cation  n-n  stacking  structure  is  also  very 
interesting.  Geometry  optimization  of  the  ionic  structure 
using  Hartree-Fock,  which  does  not  contain  electron 
correlation,  results  in  an  “open-book”  structure  in  which 
the  two  1,2,4-triazolium  planes  form  a  dihedral  angle  of 
-120°  (no  tz-tz  stacking).  Clearly,  electron  correlation 
energy  (a  major  part  of  the  aromatic  tz-tz  attraction)  in  the 
MP2  method  is  the  driving  force  for  the  formation  of  the 
cation-cation  tz-tz  stacking  structure  from  the  open-book 
non-stacking  structure. 

Figure  1.  MP2/aug-cc-pVDZ  optimized  structures  of 
two  pairs  of  1,2,4-triazolium  (1,2,4-triazole)  and 
dinitramide  (dinitramine)  molecules.  H  is  white,  C  is 
gray,  O  is  red,  N  is  blue. 


density  functional  theory  (B3LYP/6-311++G(d,p)) 
calculations. [11]  The  MP2  predicted  structures  generally 
favor  interactions  between  the  anion  and  the  n  electron 
density  of  the  cation  ring,  whereas  the  DFT  geometries 
tend  to  favor  interionic  hydrogen  bonding  between  the 
hydrogen  atoms  on  the  cation  and  the  nitrogen  or  oxygen 
atoms  of  the  anion.  An  example  of  this  difference  in 
predicted  structures  is  illustrated  in  Figure  2,  which  shows 
the  MP2  and  DFT  geometries  of  [Emim+][Im']2.  Panels 
(a)  and  (b)  ((c)  and  (d))  show  two  different  views  of  the 
same  B3FYP  (MP2)  local  minimum.  Panels  (a)  and  (b) 
show  a  predominance  of  hydrogen  bonding  interactions 
between  the  hydrogen  atoms  on  the  [Emim+]  cation  and 
the  oxygen  atoms  on  the  [Im  ]  anions.  In  contrast,  the 
MP2  optimized  geometry  of  [Emim+][Im']2  shown  in 
panels  (c)  and  (d)  suggests  the  presence  of  interactions 
between  the  anions  and  the  tz  electron  density  of  the 
cation  ring. 

Figure  2.  Calculated  geometries  of  [Emim+][Im']2, 
B3LYP/6-311++G(d,p)  (panels  (a)  and  (b))  and  MP2/6- 
311++G(d,p)  (panels  (c)  and  (d)).  H,  C,  N,  O,  F,  and  S 
atoms  are  shown  in  white,  black,  blue,  red,  green,  and 
yellow,  respectively. 
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rEmim+irim'l  Ion  Clusters:  Ionic  liquids  such  as  1- 
ethyl-3  -methylimidazolium  /?A(trifluoromethylsulfonyl) 
imide  ([Emim+][Im  ])  show  great  promise  as  advanced 
working  fluids  in  electric  propulsion  applications.  In 
order  to  assess  and  tune  the  performance  of  these  types  of 
working  fluids,  it  is  necessary  to  characterize  the  exhaust 
plume,  which  contains  a  distribution  of  clusters  such  as 
[Emim+]m[Im']m±i;  i.e.,  containing  either  an  excess  cation 
or  anion.  The  gas-phase  structures,  harmonic  vibrational 
frequencies,  and  binding  energies  of  a  series  of  ion 
clusters  [Emim+]m[Im']m±i  (m=l-3)  were  computed  at  the 
MP2/6-311++G(d,p)[10]  level  and  compared  to  earlier 


Molecular  Dynamics  Simulations  of  the  Energetic 

Ionic  Liquid  l-hydroxyethyl-4-amino-l,  2,  4- 

triazolium  Nitrate  (HEATN):  Molecular  dynamics 
(MD)  simulations  have  been  performed  to  investigate  the 
structure  and  dynamics  of  an  energetic  ionic  liquid,  1- 
hydroxyethyl-4-amino- 1 ,2,4-triazolium  nitrate 

(HEATN). [12]  The  generalized  Amber  force  field  (GAFF) 
was  used,  and  an  electronically  polarizable  model  was 
further  developed  in  the  spirit  of  our  previous  work.[13]  In 
the  process  of  simulated  annealing  from  a  liquid  state  at 
475  K  down  to  a  glassy  state  at  175  K,  the  MD 
simulations  identify  a  glass -transition  temperature  region 
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at  around  250-275  K,  in  agreement  with  experiment. [14]  A 
complex  hydrogen  bond  network  was  revealed  with  the 
calculation  of  partial  radial  distribution  functions.  When 
compared  to  the  similarly  sized  1 -ethyl-3 -methyl- 1,3- 
imidazolium  nitrate  ionic  liquid,  EMIM+/N03",  a 
hydrogen  bond  network  directly  resulting  in  the  poorer 
packing  efficiency  of  ions  is  observed,  which  is 
responsible  for  the  lower  melting/glass-transition  point. 
The  structural  properties  of  the  liquid/vacuum  interface 
shows  that  there  is  vanishing  layering  at  the  interface,  in 
accordance  with  the  poor  ion  packing.  The  effects  of 
electronic  polarization  on  the  self-diffusion,  viscosity,  and 
surface  tension  of  HEATN  are  found  to  be  significant,  in 
agreement  with  an  earlier  study  on  EMIM+/N03".[15] 

Based  on  a  charge-removed  force-matching  procedure, 
a  five-site  coarse  graining  (CG)  model  has  been 
constructed  for  HEATN.  This  model  produces  structural 
properties  which  are  in  good  agreement  with  atomistic 
simulations  (see  Figure  3).  The  coarse  graining  results  are 
also  consistent  with  the  strong  screening  effect  in  bulk 
HEATN.  A  striking  feature  of  this  CG  model  is  spatial 
transferability.  The  force  field  developed  by  a  system  of 
216  ion  pairs  has  been  successfully  used  in  larger  systems 
(343  and  512  pairs)  to  provide  good  structure  properties 
as  shown  in  Figure  3.[16]  Another  important  point  is  the 
sharp  peak  of  cation-anion  correlation  resulting  from  the 
CG  simulation:  the  double  peaks  from  the  atomic 
simulation  disappear  because  the  coarse-graining 
procedure  misses  the  oriented  information  involved  in  the 
highly  delocalized  -NH2####0N02’  and  -0H####0N02" 
hydrogen  bonding  interactions. [12,16] 


Figure  3.  Radial  distribution  functions  g(r)  of 
HEATN  from  molecular  dynamics  simulations 
utilizing  atomistic  and  coarse-graining  force  field 
models. 


216  ion  pairs 


343  ion  pairs 


512  ion  pairs 


4.  Summary  and  Conclusions 

The  geometries  of  the  dimers  of  1,2,4-triazolium 
dinitramide  and  1,2,4-triazolium  chloride  have  been 
computed  at  the  MP2/aug-cc-pvdz  level  of  theory.  These 
calculations  suggest  that  cation-cation  n-n  stacking 
structures  can  exist  in  very  small  ionic  clusters.  Electron 
correlation  energy  in  the  MP2  method  is  the  driving  force 
for  the  formation  of  the  cation-cation  n-n  stacking 
structure  from  the  open-book  non-stacking  structure.  The 
structure  motifs  and  interaction  patterns  provide  new 
understanding  of  ionic  materials  with  aromatic  rings. 

The  structures,  harmonic  vibrational  frequencies, 
binding  energies,  and  internal  thermal  energy 
distributions  of  gas-phase  ion  clusters  of  the  l-ethyl-3- 
methylimidazolium  /?A(trifluoromethylsulfonyl)imide 
ionic  liquid  have  been  predicted  using  second  order 
perturbation  theory  (MP2/6-311++G(d,p))  and  compared 
to  prior  density  functional  theory  (B3LYP/6- 
311++G(d,p))  results.  The  MP2  predicted  structures 
generally  favor  interactions  between  the  anion  and  the  n 
electron  density  of  the  cation  ring,  whereas  the  DFT 
geometries  tend  to  favor  interionic  hydrogen  bonding 
interactions. 

Classical  molecular  dynamics  simulations  of  HEATN 
predict  the  presence  of  a  glass  transition  at  a  temperature 
of  about  275  K,  in  agreement  with  experiment. 
Utilization  of  a  polarizable  versus  a  nonpolarizable  force 
field  resulted  in  significant  differences  in  predicted 
transport  properties  such  as  self-diffusion  coefficients  and 
shear  viscosity.  Simulations  based  upon  a  five-site  CG 
model  predict  structural  properties  in  good  agreement 
with  atomistic  results  and  also  show  spatial  transferability 
of  the  CG  model  to  larger  systems. 
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